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We compute the 0{o? s ) correction to the decay rate n c — > 77 and discuss its implications for 
precision quarkonium physics. We study the suitability of the ratio V(J/ip — > e + e~)/r(r;c — > 77), in 
which the non-perturbative or soft effects cancel at (a,' 1 ) , for extracting fundamental parameters 
of QCD at low energies. We show that the QCD-based theory of charmonia is not capable of 
predicting this ratio with any degree of confidence. 
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The physics of cc mesons (charmonium) is a mature 
field with a long history. Discovery of the J/ ip resonance 
at SLAC and Brookhaven in the autumn of 1974 is often 
called "November revolution" , to emphasize its impor- 
tance for the development of QCD and of the Standard 
Model. Already the early theoretical papers on the sub- 
ject interpreted the observed narrow resonance as the 
non-relativistic bound state of cc quarks and thus initi- 
ated studies of heavy quarks in QCD JjJ . 

During almost thirty years since the discovery of the 
first charmonium, both experimental and theoretical 
studies of these mesons have been successfully pushed 
forward. The spectrum, lifetimes and branching ratios 
have been precisely measured. Further progress is ex- 
pected at the planned dedicated facility CLEO-c . On 
the theoretical side, various attempts have been made to 
improve the description of these hadrons. In particular a 
lot of effort went into determining how well the cc bound 
states can be described if one starts directly from the 
QCD Lagrangian. 

An important recent development has been the in- 
troduction of effective field theory techniques for de- 
scribing hadrons consisting of two non-relativistic heavy 
quarks || Q]. This effective field theory, the Non- 
relativistic Quantum Chromodynamics (NRQCD), con- 
nects the original QCD and the new NRQCD Lagrangian 
which takes full advantage of the fact that the quarks in 
the cc bound state are non-relativistic. As usual in ef- 
fective field theories, the two Lagrangians are matched 
perturbatively at energy scales around the charm quark 
mass m c = m ~ 1.7 GeV. 

In recent years the effective field theory approach to 
non-relativistic bound states has been extended further. 



It has been noticed that two additional scales, the heavy 
quark momentum mv and the heavy quark binding en- 
ergy mv 2 , exist in quarkonia and, for sufficiently heavy 
quarks, permit perturbative treatment. Unfortunately, 
this is not quite possible for charmonium because the 
c quark mass is too small and therefore mv ~ mv 2 ~ 
Aqcd- However, we still have Aqcd "C m and there- 
fore there is a chance that integrating out hard (k ~ m) 
modes and matching QCD at NRQCD perturbatively is 
a sensible thing to do. If the soft effects are universal, 
they cancel in ratios of various observables and a clean 
perturbative QCD prediction emerges. 

For various observables, this approach has been taken 
at order 0(a s ) and the common perception is that it 
works rather well. Let us consider the simplest decays of 
the ground state charmonia, J/tp — > e + e~ and r/ c — > 77. 
To order 0(a s ) the decay rates can be written as 



il ~i/) 2 {r = 0)(l + x i ,-—} , 
,~^ 2 (r = 0) (l+av^), 



(1) 



where tp(r) is the charmonium wave function and 
are numbers which can be determined by perturbative 
matching of the QCD and NRQCD Lagrangians. Taking 
the ratio of these two decay widths, one obtains a predic- 
tion that is free from any non-perturbative uncertainties, 
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and can be either compared to the data provided a s (m) 
is known or used to extract the value of a s . For example, 
the CLEO collaboration has recently || determined the 
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strong coupling constant a s (rn c ) from the branching ratio 
of r\ c — > 77 assuming that the r\ c total decay width is 
saturated by the two-gluon channel. Since the ratio (||) 
is independent of the wave functions, its comparison with 
measurements has often been considered a solid test of 
perturbative QCD. 

In the past few years the development of techniques 
for non-relativistic effective field theories has undergone 
an important transition and we can now study the next 
order in the strong coupling constant expansion. Inter- 
estingly, in doing so one encounters new conceptual dif- 
ficulties. The most important problem is that the hard 
Wilson coefficient of the operator responsible for the de- 
cay in the leading order becomes infrared divergent at 
two loops, which implies that the wave function at the 
origin becomes scale dependent. If this scale dependence 
were the same for spin triplet J/ip and singlet r/ c states, 
this would not pose a difficulty since it would cancel in 
the ratio. However, this is not the case and the diver- 
gent parts of the Wilson coefficients are spin-dependent. 
This immediately implies that with 0(a 2 ) accuracy the 
wave functions at the origin of J/ip and r\ c become dif- 
ferent and therefore the ratio of the corresponding de- 
cay widths is sensitive to some soft-scale effects. Since 
it is rather difficult to compute these effects accurately, 
the QCD-based prediction for the ratio of decay widths 
T(J/ip — > e + e~)/T(rj c — > 77) becomes much less precise 
if one computes higher order corrections in 0{a s ) — a 
somewhat paradoxical situation. This is the principal 
message we would like to get across in this Letter. 

The remaining part of this Letter is organized as fol- 
lows. We first consider hard renormalization factors of 
the non-relativistic operators responsible for the decays 
f] c —* 77 and J/ip — > e + e~. We then show how the soft 
effects are taken into account in our calculation and de- 
rive our final result for the ratio of the decay rates of 
J/ip — > e + e~ and r\ c — > 77. 

The hard renormalization factor for the spin-singlet 
decay operator has been evaluated very recently B, ex- 
tending an earlier QED result obtained for the para- 
positronium decay |M, |f| . 
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FIG. 1: Diagrams not considered in this paper: "light-by- 
light" scattering contributions. 

We did not include the diagrams shown in Fig. [l] in 
our final result. We have checked with a rough approx- 
imation that this finite and gauge-invariant subset con- 
tributes only insignificantly. 

Our result for the hard renormalization factor for the 
singlet decay operator is (we use dimensional regulariza- 



tion with D = 4 
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where a s — a s (m) and 

s 2 {n) = C f sa + Casna + N L T R s L + N H T R s H , 



SL = 36" 
s H = 0.22. 
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In the above equation \x is the factorization scale that 
separates relativistic and non-relativistic momenta in the 
NRQCD framework. 

The matching coefficient for the vector current, rele- 
vant for the decay J/ip — * e + e~ , can be found in Eqs. (10- 
15) of Ref. § (see also fo|). Diver gences in the matching 
coefficients of the two currents are different. In the MS 
renormalization scheme the ratio of the decay widths is 
(we use C F = ±C a = 3, T R = ± N L = 3, N H = 1) 
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Obviously, the wave functions at the origin become renor- 
malization scheme and factorization scale dependent and 
eq. (||) indicates that they must be different for J / ip and 
r\ c . Since they can no longer be eliminated in the ratio 
of the decay widths, one looses predictive power. Note 
also that eq. (^) involves the MS renormalized wave func- 
tions, which cannot be directly determined even on the 
lattice. Let us stress that to reach these conclusions one 
only has to integrate out relativistic degrees of freedom. 
No assumption about the dynamics of the bound state 
(besides its quantum numbers) is necessary. 

To cancel the divergences of the Wilson coefficients we 
have to calculate certain soft effects. For charmonium, 
which is not a Coulombic bound state, a completely "hon- 
est" calculation is not possible. However, we will show 
that with plausible assumptions one can determine the 
soft contribution to the ratio of two decay widths using 
experimental data on e + e~ — > cc. Two types of relativis- 
tic corrections have to be considered: to the amplitude 
and to the wave function at the origin. It turns out that 
the former are the same for J/ip — > e + e~ and r\ c — > 77 
and cancel in the ratio (for this reason they were not in- 
cluded in eq. (^J)). As for the latter, only spin-dependent 
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effects can survive in the ratio. We assume that the only 
operator in the non-relativistic Hamiltonian responsible 
for the hyperfine splitting is 
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where the Pauli matrices a and £ act respectively on 
spins of the charm quark and of the antiquark. This 
operator follows from a one-gluon exchange diagram and 
its QED analog is the hyperfine splitting operator in the 
Breit Hamiltonian. With this operator we compute the 
wave functions at the origin, 



</>|(0) = l + ...+ 
^(0) = 1 + . 
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where G(0, 0) is the reduced Green's function of the non- 
relativistic cc state, computed at the ground state energy 
E\ (at this level we neglect the difference between J/tp 
and 7/ c masses), 



G(0,0) = £ 



l^(Q)l 2 ; 

Ei — E n 



(8) 



and the prime means that the sum does not include the 
ground state J/ip (n = 1)- The dots in eq. (|?]) indicate 
that spin-independent corrections are also present but 
they drop out in the ratio of the decay widths. 

The Green's function at the origin is divergent and it 
is precisely the divergence needed to cancel that in the 
ratio of the two hard Wilson coefficients. To illustrate 
this, we first consider an academic example of ultra-heavy 
quarkonia, where binding effects can be computed in the 
Coulomb approximation. In this case, the Green's func- 
tion at the origin G(0, 0) can be extracted from Ref. |ll[] 
and reads 

G(0, 0) = -9*2** ( 1 + 3 _ l n . ( 9 ) 

47r \4e 2 fi J w 

We now substitute this result in to eq. (^) and obtain 
the ratio of the decay widths of the spin-triplet and spin- 
singlet ultra-heavy quarkonia, 
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where a s — a s {mQ). The coefficient of the a 2 lna s term, 
4Gf73 ~ 2.37, agrees with Jl|. 

We now turn to the more difficult case of charmonium, 
where one can not use the Coulomb approximation for 
the low energy dynamics. We will still use the hyperfine 
splitting operator, eq. (^^to describe the hyperfine split- 
ting; for this reason eqs. (0) are still valid. The challenge 
is to compute the Green's function G(0, 0) in dimensional 



regularization without invoking the Coulomb approxima- 
tion. 

We will perform such computation in two different 
ways dealing either with the observed cross section 
e + e~ — > cc or with a simple potential model constructed 
to describe charmonium. It will be clear that both calcu- 
lations could be improved. Here, our primary goals are 
to demonstrate how they can be carried out using simple 
approximations and to give reasonable estimates of the 
second order corrections to the ratio of the two decay 
widths. 

We first describe the calculation which utilizes exper- 
imental data on e + e~ — * cc. To this end, we separate 
the Green's function into a contribution of charmonium 
resonances and that of the continuum, 



G(0,0) = G rcs (0,0) + G^*(0,0). 



(11) 



The former is finite and can be computed using available 
data on the spectrum and e + e _ decay widths of the spin 
1 resonances, 



G res (0,0) 



l6Q 2 c ira 2 ^ M x - M n 
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-0.073 GeV 2 , 



(12) 



where we have employed the mass and width information 
on the first six ip resonances |l3| ]. 

The continuum contribution is, on the other hand, di- 
vergent. To determine this divergence using as little in- 
put information as possible, we proceed in the following 
way. We consider continuum contribution as a function 
oiE, 



G^ nt (0,0)= J2 



E n >0 
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and take the derivative with respect to E. We then solve 
the resulting differential equation and obtain 



G&* (0,0) = G c £ nt + / d£— G c £ nt (0, 0) . (14) 

dE 



We can further use the relation between R c = <r(e + e — > 
cc)/cr(e + e _ — > and the Green's functions G(0,0) 

to relate eq. ([fi]) to experimental data. We obtain 



G|f (0, 0) = G£ nt -^zJdE K(E) R C (E), (15) 

o 

where the function K{E) is given by 

Ei — Ei 



K(E) 



(16) 



(E — Ei){E — Ei) ' 
The divergence in G(0, 0) now resides in G c £ nt and the 



integral in eq. (15) is finite. For the initial condition G c £ l 
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we can choose a "deep Euclidean" point Ei — > -co, where 
perturbative calculations are justified and where G c £ nt 
can be determined with, in principle, arbitrary precision. 
We therefore see that the divergent part of the Green's 
function can be extracted from perturbative calculations 
and the finite part can be obtained from R c (other op- 
tions are to compute the finite part using a potential 
model for quarkonium or NRQCD on the lattice). This 
separation solves the problem in principle and provides 
a way to determine the MS charmonium wave functions 
without assuming that the bound state is Coulombic; the 
only true model-dependence remaining in this calculation 
is the form of the operator responsible for the hyperfine 
splitting, eq. (^|). 

The data on R c is not quite precise yet. Its high energy 
asymptotics (in the non-relativistic sense) is fixed since 
the dependence on the initial energy Ei in eq. ([lj]) should 
cancel. We therefore write 



Rc(E) 




R n c pt (E). 



For Rc pt (E) we choose: 

B£V*(E) = ~irC F a s 9(E Q - E), 



(17) 



(18) 



with Eq — y/~SQ — 2m and y/s Q = 4 GeV. This ansatz 
is motivated by the data on R c in where one sees 
that there is no need for a large second term in eq. fll7j ) 
below 4 GeV. On the other hand we do need this term 
at higher energies, since otherwise "perturbative" and 
"non-perturbative" expressions do not match. 
We obtain 



G c if (0,0) 
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The soft contributions to the decay width ratio are ob- 
tained by employing eqs. ( |l2| ) and ( |l9| ) in eq. (@). Before 
combining them with the hard contributions in (^]), we 
note that the BLM effects were computed in Ref. |l(| 
for the rates of J/tfj — ► e+e~ and rj c — > 77. These correc- 
tions turn out to have different signs and are enhanced 
in the ratio. For this reason, we decided to eliminate 
them by choosing different scales for the strong coupling 
constant in the Wilson coefficients for J/tp — > e + e~ and 
77 c — > 77. We then find 



J/ijj — >e+e 



+at 1.19 In 



1 

= 3Qf 

El — Eq 



1 - 1.7 as=i + 1-1 as=o 



3.66 




1.64 



rnttj a s m 2 



(20) 

where ag =0 = a s (1.95 m) and as=i = a s (0.63 to) and 
the charm mass should be expressed in GeV. The scale 



of the coupling constant in the second order correction is 
not specified; we will use a s — 0.3 for the estimates. The 
inverse power of a s in the square brackets arises because 
we have used experimental data to compute the contri- 
bution of the resonances to G(0,0) and also used the 
energy of the ground state E\ = Mj/^ — 2m to estimate 
the continuum contribution. This spoils the homogeneity 
in a s . 

Employing, for the sake of illustration, a s — 0.3, 
a s (0.63 to) = 0.35 and a s (1.95 to) = 0.26 in eq. (|2C 
one finds 



J/i>- 



1 



r 

1 no^n 



= — [1-0.32 + h(m)+0(a 3 s )] 



(21) 



where the three terms in the brackets are the tree level, 
the 0(a s ) and the 0{a 2 s ) corrections, respectively. As 
shown in Fig. the second order correction depends 
strongly on the value of the charm quark mass. It is 
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FIG. 2: The second order correction in eq. (^l]) as a function 
of the charm quark mass. 

large and very sensitive to the form of i?" pt assumed in 
eq. (^). For example, if we use the Coulomb approxi- 
mation to estimate the value of the wave function in the 
continuum, we obtain a negative result for the second or- 
der correction. Within our model for the continuum, we 
estimate 



~ = ^ [1 - 0.32 + 0.3 ± 0.1 + 0(a3)]. 

(22) 

As we mentioned earlier, another possibility to obtain 
the Green's function necessary to compute the ratio of 
two decay widths is to use either QCD on the lattice or 
potential models. Here we would like to illustrate this 
point by considering a simple potential model. This will 
allow us to check that the estimate, eq. (p2|), is reason- 
able. Although the potential model below is really sim- 
ple, the calculation can be repeated with more sophisti- 
cated potentials, provided that, at short distances, these 
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potentials match the QCD analog of the Coulomb poten- 
tial. 

To describe charmonium, we will use the Schrodinger 
equation with the potential 



V(r) 



C F a s 



br + Vo, 



(23) 



where b = 0.18 GeV 2 and Vb is adjusted to give the 
correct mass of the IS state (J/tp) f° r given values of m 
and a s . 

To compute the Green's function G(0, 0) we solve the 
Schrodinger equation following the treatment in jlTj] and 
obtain the following representation for the full (including 
the ground state) Green's function: 



G E (0,0) = lim- — 

r^Q 4-7T 



- - C F a s mhx{r) + B(E) 



, (24) 



where the function B{E) is derived from the large r limit 
of the ratio of two solutions of the Schrodinger equation 
with prescribed behavior at the origin (see |l7J ) . In order 
to obtain the Green's function in dimensional regulariza- 
tion, which is needed for our purposes, we write 



G_b(0,0) = Ge~Ge 1 



Ge 
m 

47T 



B(E 1 



Gi 



(25) 



We then take the limit E\ — > — oo, perturbatively com- 
pute the Green's function Ge x and derive 
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where c is 
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and je is the Euler constant. 

In order to compute the reduced Green's function 
G(0, 0) eq. (||), we find the first eigenvalue numerically 
and remove the pole from B{E). 

We have computed the ratio of two decay widths for 
different values of a s and the mass of the charm quark. 
We obtain 



1 



r 



— 2 [l - 0.32 + 0.20 ± 0.05 + O (a 3 s )] 



3Q 



(28) 



We see that the result of the potential model calculation 
is relatively close to the result of eq. (22). The advan- 
tage of the potential model calculation is its fair stability 



against variations in a s and m. We take the result in 
eq. (p8|) as our final estimate. 



In spite of the fact that the model leading to eq. (|2^ 
is quite simple, we believe that eq. ( ^8| ) is important in 
that it clearly shows the magnitude of second order QCD 
corrections one might expect for such observables. 

It is interesting to note that there is a strong can- 
cellation between the first and second order effects in 
eqs. ( p2| , p8| ). Neglecting all the radiative corrections 
and using T^„, e + e - = 5.26 keV, we derive r^ c _> 77 = 
7.01 keV, rather close to the central value reported by 
CLEO collaboration @ r^ 77 = [7.06 ± 0.8(stat) ± 
0.4(sys) ± 2.3(br)] keV. 

We conclude that eqs. ( pp| , pl| , p8| ) , the principal results 
of this Letter, illustrate an unexpected problem in the 
theory of heavy quarkonia at the two-loop level. In recent 
years we have learned how to integrate out relativistic 
degrees of freedom efficiently and it seemed as if we could 
improve the accuracy of our predictions. This turns out 
not to be the case. The reason is that at 0(a 2 s ) the 
soft and relativistic effects do not decouple completely, 
as it happens at 0(a s ), and therefore, in general, one 
cannot avoid non-perturbative effects by taking ratios of 
different observables. We have shown how, in principle, 
the soft contribution can be estimated using experimental 
data or potential models. 



With the QCD corrections as big as in eq. (|28|), the 
determination of a s (m c ) from charmonia decay rates, as 
e.g. in Ref. [^), does not look trustworthy, regardless of 
the fact that the numerical values of a s turn out to be 
in a theoretically sensible range. On the other hand, it is 
interesting to point out that the ratio of the decay rates 
of rj c — > 77 to r/ c — * gg, actually used in Ref. [^) for deter- 
mination of a s , is free from the soft effects we discussed 
in this Letter, since it refers to the same initial state. It 
would therefore be interesting to compute second order 
QCD corrections to this ratio since in this case the hard 
corrections alone might provide an unambiguous answer. 

Among various charmonium decays, only J/ip — > e + e~ 
and r\ c — > 77 have now been studied to 0(a 2 ). Clearly 
these are the two simplest channels since they do not 
involve any complications related to the dynamics of 
hadrons in the final state. If the understanding of even 
those simplest decays encounters such difficulties, one 
should exercise great care when extracting physical in- 
formation from more complicated charmonium decays. 
The fact that perturbative QCD appears to work well in 
the one-loop order is certainly insufficient to ensure that 
heavy quarkonia are well understood. 
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